﻿using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace SmartMathLibrary.NumericalIntegration
{
    /// <summary>
    /// This class provides a bivariate numerical integration by using the mid point method.
    /// </summary>
    [Serializable]
    public class BivariateMidPointIntegrator : AbstractBivariateIntegrator
    {
        /// <summary>
        /// Initializes a new instance of the <see cref="BivariateMidPointIntegrator"/> class.
        /// </summary>
        /// <param name="integral">The bivariate integral to integrate.</param>
        public BivariateMidPointIntegrator(BivariateIntegral integral)
            : base(integral)
        {
        }

        /// <summary>
        /// Integrates the integral in a numerical way. After running the numerical integration
        /// it should be checked if the result has a precision error.
        /// </summary>
        /// <param name="iterations">The maximum numbers of iterations.</param>
        /// <param name="precision">The minimum precision of the numerical integration.</param>
        /// <returns>The area of the specified integral.</returns>
        public override double Integrate(int iterations, double precision)
        {
            double a = Math.Max(this.integral.A1, this.integral.B1) - Math.Min(this.integral.A1, this.integral.B1);
            double b = Math.Max(this.integral.A2, this.integral.B2) - Math.Min(this.integral.A2, this.integral.B2);

            return this.Integrate((int) ((a * b) * 50), iterations, precision);
        }

        /// <summary>
        /// Integrates the integral in a numerical way. After running the numerical integration
        /// it should be checked if the result has a precision error.
        /// </summary>
        /// <param name="parts">The number of parts, in which the function should be divided at the start of the algorithm.</param>
        /// <param name="iterations">The maximum numbers of iterations.</param>
        /// <param name="precision">The minimum precision of the numerical integration.</param>
        /// <returns>The area of the specified integral.</returns>
        public override double Integrate(int parts, int iterations, double precision)
        {
            double sign1 = 1;
            double sign2 = 1;
            this.a1 = this.integral.A1;
            this.a2 = this.integral.A2;
            this.b1 = this.integral.B1;
            this.b2 = this.integral.B2;

            if (b1 < a1)
            {
                double tempuri = b1;

                b1 = a1;
                a1 = tempuri;
                sign1 = -1;
            }

            if (b2 < a2)
            {
                double tempuri = b2;

                b2 = a2;
                a2 = tempuri;
                sign2 = -1;
            }

            int i = parts;
            double newValue = 0;
            double oldValue = this.InternalIntegrate(i);

            iterations += parts;

            while (i < iterations)
            {
                i++;
                newValue = this.InternalIntegrate(i);
                this.relativeError = Math.Abs(newValue - oldValue);

                if (this.relativeError < precision)
                {
                    this.neededIterations = i;
                    this.precisionError = false;

                    if (sign1 == sign2)
                    {
                        return newValue;
                    }

                    return -newValue;
                }

                oldValue = newValue;
                //i++;
            }

            this.neededIterations = i;
            this.precisionError = true;

            if (sign1 == sign2)
            {
                return newValue;
            }

            return -newValue;
        }

        /// <summary>
        /// Internal integration algorithm for the used integration method.
        /// </summary>
        /// <param name="parts">The number of parts to divide the area in.</param>
        /// <returns>The integrated area of the specified integral.</returns>
        private double InternalIntegrate(int parts)
        {
            double tempuriSecondDimension = 0;
            double hSecondDimension = (this.b2 - this.a2) / (double) parts;

            for (int i = 0; i < parts; i++)
            {
                double tempuriFirstDimension = 0.0;
                double hFirstDimension = (this.b1 - this.a1) / (double) parts;

                for (int j = 0; j < parts; j++)
                {
                    tempuriFirstDimension += this.integral.Function.SolveAt(this.a2 + (i + 0.5) * hSecondDimension,
                                                                            this.a1 + (j + 0.5) * hFirstDimension);
                }

                tempuriSecondDimension += hFirstDimension * tempuriFirstDimension;
            }

            return hSecondDimension * tempuriSecondDimension;
        }
    }
}